42 research outputs found
Rational Krylov for Stieltjes matrix functions: convergence and pole selection
Evaluating the action of a matrix function on a vector, that is , is an ubiquitous task in applications. When is large, one
usually relies on Krylov projection methods. In this paper, we provide
effective choices for the poles of the rational Krylov method for approximating
when is either Cauchy-Stieltjes or Laplace-Stieltjes (or, which is
equivalent, completely monotonic) and is a positive definite
matrix. Relying on the same tools used to analyze the generic situation, we
then focus on the case , and
obtained vectorizing a low-rank matrix; this finds application, for instance,
in solving fractional diffusion equation on two-dimensional tensor grids. We
see how to leverage tensorized Krylov subspaces to exploit the Kronecker
structure and we introduce an error analysis for the numerical approximation of
. Pole selection strategies with explicit convergence bounds are given also
in this case
A low-rank technique for computing the quasi-stationary distribution of subcritical Galton-Watson processes
We present a new algorithm for computing the quasi-stationary distribution of
subcritical Galton--Watson branching processes. This algorithm is based on a
particular discretization of a well-known functional equation that
characterizes the quasi-stationary distribution of these processes. We provide
a theoretical analysis of the approximate low-rank structure that stems from
this discretization, and we extend the procedure to multitype branching
processes. We use numerical examples to demonstrate that our algorithm is both
more accurate and more efficient than other approaches
Solving rank structured Sylvester and Lyapunov equations
We consider the problem of efficiently solving Sylvester and Lyapunov
equations of medium and large scale, in case of rank-structured data, i.e.,
when the coefficient matrices and the right-hand side have low-rank
off-diagonal blocks. This comprises problems with banded data, recently studied
by Haber and Verhaegen in "Sparse solution of the Lyapunov equation for
large-scale interconnected systems", Automatica, 2016, and by Palitta and
Simoncini in "Numerical methods for large-scale Lyapunov equations with
symmetric banded data", SISC, 2018, which often arise in the discretization of
elliptic PDEs.
We show that, under suitable assumptions, the quasiseparable structure is
guaranteed to be numerically present in the solution, and explicit novel
estimates of the numerical rank of the off-diagonal blocks are provided.
Efficient solution schemes that rely on the technology of hierarchical
matrices are described, and several numerical experiments confirm the
applicability and efficiency of the approaches. We develop a MATLAB toolbox
that allows easy replication of the experiments and a ready-to-use interface
for the solvers. The performances of the different approaches are compared, and
we show that the new methods described are efficient on several classes of
relevant problems
Low-rank updates and a divide-and-conquer method for linear matrix equations
Linear matrix equations, such as the Sylvester and Lyapunov equations, play
an important role in various applications, including the stability analysis and
dimensionality reduction of linear dynamical control systems and the solution
of partial differential equations. In this work, we present and analyze a new
algorithm, based on tensorized Krylov subspaces, for quickly updating the
solution of such a matrix equation when its coefficients undergo low-rank
changes. We demonstrate how our algorithm can be utilized to accelerate the
Newton method for solving continuous-time algebraic Riccati equations. Our
algorithm also forms the basis of a new divide-and-conquer approach for linear
matrix equations with coefficients that feature hierarchical low-rank
structure, such as HODLR, HSS, and banded matrices. Numerical experiments
demonstrate the advantages of divide-and-conquer over existing approaches, in
terms of computational time and memory consumption
On maximum volume submatrices and cross approximation for symmetric semidefinite and diagonally dominant matrices
The problem of finding a submatrix of maximum volume of a matrix
is of interest in a variety of applications. For example, it yields a
quasi-best low-rank approximation constructed from the rows and columns of .
We show that such a submatrix can always be chosen to be a principal submatrix
if is symmetric semidefinite or diagonally dominant. Then we analyze the
low-rank approximation error returned by a greedy method for volume
maximization, cross approximation with complete pivoting. Our bound for general
matrices extends an existing result for symmetric semidefinite matrices and
yields new error estimates for diagonally dominant matrices. In particular, for
doubly diagonally dominant matrices the error is shown to remain within a
modest factor of the best approximation error. We also illustrate how the
application of our results to cross approximation for functions leads to new
and better convergence results
On Functions of quasi Toeplitz matrices
Let be a complex valued continuous
function, defined for , such that
. Consider the semi-infinite Toeplitz
matrix associated with the symbol
such that . A quasi-Toeplitz matrix associated with the
continuous symbol is a matrix of the form where
, , and is called a
CQT-matrix. Given a function and a CQT matrix , we provide conditions
under which is well defined and is a CQT matrix. Moreover, we introduce
a parametrization of CQT matrices and algorithms for the computation of .
We treat the case where is assigned in terms of power series and the
case where is defined in terms of a Cauchy integral. This analysis is
applied also to finite matrices which can be written as the sum of a Toeplitz
matrix and of a low rank correction
Efficient cyclic reduction for QBDs with rank structured blocks
We provide effective algorithms for solving block tridiagonal block Toeplitz
systems with quasiseparable blocks, as well as quadratic matrix
equations with quasiseparable coefficients, based on cyclic
reduction and on the technology of rank-structured matrices. The algorithms
rely on the exponential decay of the singular values of the off-diagonal
submatrices generated by cyclic reduction. We provide a formal proof of this
decay in the Markovian framework. The results of the numerical experiments that
we report confirm a significant speed up over the general algorithms, already
starting with the moderately small size
A nested divide-and-conquer method for tensor Sylvester equations with positive definite hierarchically semiseparable coefficients
Linear systems with a tensor product structure arise naturally when
considering the discretization of Laplace type differential equations or, more
generally, multidimensional operators with separable coefficients. In this
work, we focus on the numerical solution of linear systems of the form where the matrices are
symmetric positive definite and belong to the class of hierarchically
semiseparable matrices.
We propose and analyze a nested divide-and-conquer scheme, based on the
technology of low-rank updates, that attains the quasi-optimal computational
cost where is the condition number of the linear
system, and the target accuracy. Our theoretical analysis highlights
the role of inexactness in the nested calls of our algorithm and provides worst
case estimates for the amplification of the residual norm. The performances are
validated on 2D and 3D case studies
Optimizing network robustness via Krylov subspaces
We consider the problem of attaining either the maximal increase or reduction
of the robustness of a complex network by means of a bounded modification of a
subset of the edge weights. We propose two novel strategies combining Krylov
subspace approximations with a greedy scheme and an interior point method
employing either the Hessian or its approximation computed via the
limited-memory Broyden-Fletcher-Goldfarb-Shanno algorithm (L-BFGS). The paper
discusses the computational and modeling aspects of our methodology and
illustrates the various optimization problems on networks that can be addressed
within the proposed framework. Finally, in the numerical experiments we compare
the performances of our algorithms with state-of-the-art techniques on
synthetic and real-world networks